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Multidisciplinary Software 


I. Introduction 

The Thaerocomp Technical Corporation (TTC) has proposed to develop next-generation software 
that will revolutionize practical CFD-based analyses by providing the capability for seamless 
multidisciplinary simulation using highly accurate procedures. The software identifies, 
integrates, and then systematically upgrades common elements of major research codes, to yield a 
versatile high-order (Fele, 1991) software suite that can be used to simulate, in a coupled fashion, 
various elements of turbulence, electromagnetics, aeroacoustics, magnetogasdynamics, fluid- 
structure interaction, and chemically reacting flows. An elegant and intuitive graphical user 
interface (GUI) will be developed by leveraging experience derived from our flagship product 
INSTED. TTC’s new, tightly integrated and versatile product will reduce the simulation time by 
orders of magnitude and drastically reduce the design cost of modern engineering systems 
requiring multidisciplinary optimization. The proposed software will facilitate out-of-the-box 
approaches by providing a validated realistic virtual environment, thus impacting all major thrusts 
of the Air Vehicles Directorate, including access-to-space, sustainment, and unmanned air 
vehicles. 

The present document summarizes the accomplishments of the Phase I effort, which successfully 
and firmly established the feasibility of the concept proposed by TTC. In specifics, the following 
tasks were successfully carried out in Phase I: 

• Identification of several existing, but not necessarily high order, upgradable target codes, with 
a synopsis of the numerical characteristics of the codes. For each code, several phone and 
written contacts were made in order to obtain the current status of the code. 

• Selection of the base code for the development of the high order multidisciplinary software. 
The specific criteria for selection include a) the type of numerical method (finite 
element/volume, finite difference, spectral methods), b) capabilities of the codes, c) grid 
topology (structured/unstructured), d) the structure of the code, vis-a-vis the incorporation of 
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chimera/overset capabilities for competitive modeling of complex geometry and relative 
motion, e) the availability of the code. 

• Identification of the numerical kernels in the codes and development of a generalized flux 
function that supports all the applications. 

• Implementation of an elegant and intuitive graphics user interface (GUI) to demonstrate the 
feasibility of a user-friendly interface for the high order multidisciplinary software, which is 
an absolute necessity for commercialization of the final product. 

• Development and testing of a preliminary version of the proposed software. A demonstration 
for MGD and Aeroacoustics is presented in terms of the computational efficiency and 
accuracy of computation when compared with codes that are dedicated to these applications. 

• Investigation of marketing efforts for the commercial sector during Phase I, which shows that 
considerable interest exists for the proposed product because of its potential to ensure fast 
time-to-market. 

II: Identification of Target Codes 

The target codes for upgrade have been identified, as were the kernels in these codes. Some of 
these codes are listed in the table below, together with a synopsis of their numerical 
characteristics. Several contacts were made to find out the current state of the various codes. Dr. 
Pieter Buning (NASA Langley) provided an overview, current status, and the relationship 
between ARC3D code (originally developed by Tom Pulliam), XAIR (derivative of ARC3D), 
OVERFLOW (also a derivative of ARC3D), CFL3D (developed from NASA Langley), and 
WIND (originally developed by Douglas Ray Cosner/McDonnell Douglas). Greg Power 
(Contractor at AEDC) provided detailed information to TTC on the current status of WIND, 
which is the main production flow solver for the NPARC Alliance. 



ARC3D 

FDL3DI 

XAIR 

OVERFLOW 

WIND (NPARC) 


(NASA) 

(WPAFB) 

(AEDC) 

(NASA) 

(Boeing) 







Method 

FD 

FD 

FD 

FD 

FD (Node-based 

FV) 

Spatial 

Scheme 

Upwind, 

centered 

Compact, 

centered, 

upwind 

Upwind, 

centered 

Upwind, centered 

Upwind, centered, 

Spatial 

accuracy 

2 nd , probably 

3 rd (Roe) 

2 ,4,6, 
etc. 

2 nd , probably 
3 rd (Roe) 

2 nd , probably 3 rd 
(Roe) 

2 nd , probably 3 rd 
(Roe) 

Temporal 

scheme 

Beam- 

Warming 

Beam- 

Warming, 

RK4 

Partial 
factorization 
+ Gauss- 
Seidel 

Beam-Warming 

Partial factorization 
+ point Jacobi or 
Gauss-Seidel 

Time 

Accuracy 

Steady state, 
second order 
with sub¬ 
iterations 

RK4, 2 nd , 

3 rd , time - 
accurate 

2 nd , time- 
accurate 
moving 
grids 

Steady state, 

2 nd -order with 
sub-iterations 

First-order, 

essentially 

Stabilization 

Explicit 4 th /2 nd 
dissipation, 
Inherent Roe 

Compact 
filters, 
explicit 4/2 

HLLE (Roe¬ 
like) 

upwinding 

Explicit 4 th /2 nd 
dissipation, 
inherent Roe 

Upwinding 

Multiblock 

No 

Yes 

Yes 

Yes 

Yes 

Chimera, 

Overset 

Capabilities 

None 

Limited 

Extensive 

Extensive 

Extensive 
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Turbulence 

modeling 

RANS 

LES 

RANS 

RANS 

RANS 

Matrix 

method 

Scalarized 

pentadiagonal 

TDMA 

N/A 

Scalarized 

pentadiagonal 

N/A 

Mesh 

interpolation 

tools 

Limited 

Limited 

Extensive 

Extensive 

Extensive 

Extent of 
dissipation 

Significant 

Slight 

Significant 

Significant 

Significant 

Applications 

Standard 

Navier-Stokes, 
research, 
mostly in- 
house 

Acoustics, 
DNS/LES, 
CEM, MHD, 
mostly in- 
house 

In-house 
only, not 
supported 

Extensive industry 
use 

Part of NPARC, 
extensive industry 
use 


On the basis of the numerical characteristics of the reviewed code, vis-a-vis the goal of this 
project, we selected the FDL3D code (Gaitonde and Visbal, 1998) as the basis for the 
development of the multidisciplinary CFD tool. The specific criteria used include the type of 
numerical method, the capabilities of the target code, grid topology in the code, chimera/overset 
capability and the availability of the code. Note that none of the versions of the FDL3D code is 
capable of multidisciplinary simulation. They also do not currently provide high-order shock¬ 
capturing capabilities. Furthermore, the overset/chimera capability is not as advanced as that in 
some of the other codes that we reviewed, notably OVERFLOW and WIND. The types of 
boundary conditions currently supported in FDL3D are limited, for example, when compared 
with WIND. These 4 areas were automatically added to the technical focus for the SBIR project. 
The WENO scheme (Shu, 1997) was selected for high order shock capturing. 

Ill: Governing Equations and mathematical model 

The Phase I proposal promised to develop a preliminary generalized framework for the 
multidisciplinary software, including a generalized flux function, and to select the platform for 
the high order calculations and for stabilization in flow fields with shocks and/or strong 
discontinuity. It was also planned to test a preliminary version of the software suite. 

A preliminary version of the proposed module that supports generalized flux functions and source 
terms was developed in Phase I. The module uses the same high order procedure for evaluating 
the derivatives for the various flux functions (source terms). The analysis of the efficiency and 
accuracy of the generalized flux procedure relative that of the individual single applications was 
also started in Phase I. Each of the supported applications, in curvilinear coordinates, can be 
written in the form 


8X f dF 8G f 8H f 
dr + + dr] + d£ 



where x is the time domain, (% f T), Q are the transformed coordinates that respectively correspond 
to the physical coordinates (x,y,z). F\ G\ and H’ are the inviscid fluxes, while F’ v , G’ v , and H 9 v 
are the fluxes associated with the stresses for a particular application. S’ is a generic source term. 
Some manipulations are sometimes required in order convert the stresses into the foregoing form, 
notable examples being the Maxwell's stresses in MGD. The various flux functions that need to 
be “combined” are listed below, using F” as an example. The generalized flux functions for some 
of the multidisciplinary applications are presented subsequently. 
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DNS , non-MHD : When a handle to the module indicates direct numerical simulation , or the fluid 
component of a fluid-structure interaction , the module assumes the following definitions 



p 


( pU \ 

pu 


puU + £ x p 

pv 

F = 

pvU + £ y p 

pw 


pwU + ^p 

pE 


\ p(E + p)U —£tP J 


where 


= U +5^ + ^ +4*u 

LES, non-MHD : For non-MHD, LES problems, the following definitions apply (Germano et ah, 
1991; Rizetta et al, 2000, 2001): 


9 


( ~ ptJ \ 

pu 


puU + ZsP 

-pu 

r = 

pvU' + £ il p 

pw 

."lW 


pwU + 

pE 


\ pEU + UiP J 


where 




MGD \ For a simplified MGD model, i. e., assuming chemical and thermal equilibrium, the 
appropriate definitions are (Gaitonde, 1999): 




pU \ 

pvU + £j> -£ x RiB a 2*- 

r v r 'y,\ 

P»U + { a p-t.R i B„fL 
ptut7 + 0> 

_ r'Jvi rwi r'-'Ki _ 

(pE+p)U -^^(U • f-)B. -iMJJ- f-)B s 

r'y,\ " r v r,\ r'-'Ki 

^ - uB s ) 

iJuB y -vB^ + ^jwBy -vB s ) 

£J l uB s -wB^ + ZytyBn -wB y ) j 


The standard notations apply (see reference). 

Species conservation: When reacting or non-reacting species have to be simulated, the module 
must branch to the following definitions: 



’ 9<t> i " 


( 9<PiU \ 


9<P 2 

, F* = 

p<p 2 U 


. P#* . 


V P^ U / 
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where, as before, we can define 


u = £* + + + 

Note that the species fluxes could easily be appended to any of the foregoing definitions, in a 
coupled approach, as opposed to the segregated manner implied here. 

CEM: For the time-dependent electromagnetics, the modules sends controls to the Maxwell's 
equations for the electromagnetic field in free-space, which can be written as (Shang and 
Gaitonde, 1995, Gaitonde and Shang, 1997): 



0 


{ ^D s /e-^Dje \ 


0 


-inDz/e + ^Dv/e 

1 

0 

P f — 

^Dy/e-^Dv/e 

J 

J, 

1 r ~ 



Jy 


^n;B z If^m //Crt 


J , 


\ -ZxBy/Vn+ZvBx/Vn J 


In these expressions, B = (B x , B y , B z ) is the vector of magnetic field density, D = (D x , D y , D z ) is 
the vector of electric displacement, p m is magnetic permeability, and e is the electric permittivity. 

The starting point for the completion of the generalized flux function in Phase II is the assembled 
functions, which are defined below: 


dx f dF dG f dH f 

dr + + drj + dC; 


= (1 — 7T 6 -7T T ) 


-dFi dGi dHj- 
drj <9C 


+7T2S2 + = S, 



DNS 

LES 

CRF 

MHD 

CEM 

CAA 

FSI 

71 Fi 

1 

0 

0 

0 

0 

0 

0 

rt 2,j 

0 

1 

0 

0 

0 

0 

0 

rt 3,j 

0 

0 

1 

0 

0 

0 

0 

rc 4,i 

0 

0 

0 

1 

0 

0 

0 


0 

0 

0 

1 

0 

0 

0 

rt 6,i 

0 

0 

0 

0 

1 

0 

0 

rt 7,j 

0 

0 

0 

0 

0 

1 

0 

JE8J_ 

0 

0 

0 

0 

0 

0 

1 


In this table, the subscript “j” on n has been introduced to identify the columns. The abbreviations 
in the table are as follows: 


Symbol 

Description 

Remarks 

DNS 

Direct numerical simulation 

Standard Navier-Stokes Eqns 

LES 

Large eddy simulation 

Standard Navier-Stokes Eqns 

CRF 

Chemically-reacting flow 

Continuum 

MGD 

Magnetogasdynamics 

Continuum 

CEM 

Computational electromagnetics 
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CAA 

Computational aeroacoustics 


FSI 

Fluid-structure interaction 

Structure equations will be 
included in the Phase II 
project 


The complete flux functions for each of the applications are subsets of the following generalized 
functions derived during Phase I: 


x “ j 


(tTi + -7T2 + -TTg + 7T4 + rr E + 7T T}p + TteBn 
(tTi + 7T2 + TTs + 7T4 + TTs + TTj )pU + W^By 

(tti + tt 2 + TTs + 7T4 + rr E + 7 Tt)^u + 

(tTi + 7T2 + -TTg + 7T4 + TTe + + 7r 6t 7 K 

(tTi + 7 T2 + TTg + 7Tr )p + (^4 + TTe )pEs + W^Jy 

■TTg pYi + 7T 6 D S + 7T e J s 
(^4 + ^)5^ 

(rr 4 + 7t b )B u 
(tt 4 + tt b )B s 


The inviscid fluxes are 


/ (-7T1 +7T2 +7Tg + 7T 4 + TTfi + 7T T )pU + 7Te(^B s /e - ^D^/e) \ 

K + 7T 2 + TTg + 7T 7 ) \puU + + (tT 4 + 7T S ) - R^f*-] 

r ■ s-i-i 

+7r 6 (-^D :£ /e + 4D K /e) 

(TTx + 7T 2 + 7T 3 + 7T T ) [0UI7 + £ p] + (tT 4 + 7T E )B 

r'-' k-i 

+7r 6 (^D y /e -^Dx/e) 

(■TTl +7T2 +7Tg + 7T T )[|WU 17 + 4p] + (^4 + TTfi)[4 E b - ^5^] 

+ 4^/0 

(■7T]_ + ^2 + ^3 + - 717 ) (pi? +p)f7 + (tt 4 + tt e ) [2i?;,t7 — R&(U ■ B)^ 0 -] 

1 "' m 

Trg^yilT] +^5,/^) 

(tt 4 + tt £ )+ w^)B x - u{£ y B y + 45*)] 

\ (-7T4 + TTe) [(^4 + V^ y )B s — wl^Bn + / 


and 


/ 



V 


(TTi +7T2 +7Tg +-7T4 + -7rE +7T 7 )pV + 7T 6 D s /e - 7} s Dy /V) 

fri + rr 2 + rr 3 + tt 7 ) [/mV + tjj?] + (tt 4 + 7r e ) - R^f^] 

r' m 

+7r 6 (-^Z) s /e + 7] :£ r) 5= /e) 

K + 7T 2 + 7T 3 + 7T T ) [puV + 7^/?] + (tT 4 + TTe)^^ - R b B yf^] 

+MVvDyi£-'%Dv/£) 

(?Ti + 7T2 + TTg + 7T T )[^ V + 7^/?] + (tT 4 + 7Te)[^X & - RjBy 

i" m 

+M~Vy B */p’ m +V* B yfa m ) 

(tTi + 7T2 + TTg + 7T T )(p-E + ^)V + (tT 4 + TTfi) [2EbV ~ R&(U ■ B)f“-] 

r' ki 

+7T6 (ij^a //* m - //O 

TrabliV] +^ 6 (-vA/^ +'%Bn/^) 

(tt 4 + 7T fi )[u^ + wri^Bx - u(rj y By +i} x B x )] 

(^4 + ^s)[(^ + -ufaA + T^R*)] 

(7r 4 + 7r s )[(u77 aj + u77 s )R s F^R*)] 


\ 
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/ (TTi +7T2 +7T3 +7T4 + 7TS + ^ T )pW + n 6 (CyB s /e - QDJe) \ 

(TTi +7T 2 + 7T 3 +7T T )[^iy + C^] + (tT 4 + 7T e ) [C E ~ 

r- s-i-i 

( — C As/e + C s As /e) 

(tTi + 7T 2 + -TTg + 7T 7 ) \pvW + Q r p] + (tT 4 + 7T S ) [CyB b ~ RbB s j^-] 

■* r' 1 ' r.\ 

+^{Q s D y ie -Q y D s ie) 

frl + 7T 2 + 7T 3 + n Y )[pwW + Cjj] + (7T 4 + 7T £ ) [C s ^ “ 
+'K6(-£ y B s fp, m + Cs^/aO 

(TTi + 7T 2 + 7T 3 + 7T T )(^^ + p)W + (rr 4 + K E )\ 2 E h W - jRj, (U ■ B)f*-] 

ip>m — 

Tr^riiy] + 7T 6 (-CB y /^ m + c y ^ K /^ m ) 

(^4 + ^)^^ +^C)£ k - u(CjA + CA)] 

(^4 + TTe) [« e + wtJBy - u (C e B k + CA)\ 

\ (n 4 + n E )[(uC x +vCy)B^-wiC^B^ + CyBy)] } 


The dissipative fluxes are as follows: 


( 


Ft. 


1 

J 


0 \ 

[frl + ^2 + ^3 + ^4 + 7T S )Tii + (7T2 + 

4-Ja* [(^l + ^2 + 7T 3 + 7T 4 + 7T S )Ti2 + (^2 + 

B^a* [(^1 + K 2 + TTs + 7T 4 + 7Ts)Ti3 + (tT 2 + 

(£k, ^ [(TTi + 7T2 + 7T3 + 7T 4 + 7 T S )t^ + (tT 2 + ^ 3 )^/®] 

+ fr-l)FrM» [(^1 + 7T 2 + 7T 0 + 7T 4 + 7T B )® + (7T 2 + 7Tg)Qf £S ] 

+(rr 4 + tt e )- 


.■Fl-l, 






3 s. 


^IjO- 

1 


(tt 4 + TTs) 
(tT 4 + 


f -LfiSll -£ f 

® s -f '■2jU ffl ^ 3a; ju 

' ^7^k3(^ + %f/ S )] 

_|_ £ _3_JLil. _i _3_J?: 
~ 3s *y £>& ^ 

i_s JLAl _i _3_s, 

T W 3* u... 


i .3.^ 
^ 3y 

£ £_ElL 

“5s 3 k 

£ i.5*. 

w a® 


1 _£ Sf 


+ 


j?_ 

3s «■„ 


_3_^. 

13 Pm 

+ £ AE^ 
^>>v3yft m 


_e _3_S, 
3& fA, 


_p 

_£ d-M*. 

_i 3 ^ 

‘Sy as 


■)' 




( 


(31. = 4 


1 

j 


0 

[(^1 + ^2 + TTS + -7T4 + TTfiJm + (tT 2 + ^ 3 )^^] 

[C^I + rr 2 + 7T 3 + 7T 4 + 7r s )T i2 + (rr 2 + ^)^f 2 ES ] 

Rl^K, [(^1 + n 2 + ?r 3 + n 4 + 7T E )Ti 3 + (tt 2 + ■7r 3 )crfg E '' s 'l 

[(^1 + -7T2 + 7T 3 + 7T 4 + TTsjTij + (?T 2 + 7T 3 )fff/^] 

[(^1 +7r2 + 7r 3 + 7r 4 + 7r s )$ + (tt 2 +^)Qf ES ] 


3 Bi 


+(7r 4 +7r E ) ? 3^ ka^Ca^) 

K®« + ^/ s )l 


1 _£ Si 

1 ^ 


I 




(■7T 4 + TT S ) 1 

(tT 4 + -TTs ) R 
^4+7T E )^ ff 


jRc 'j - 
1 

tr O’ 
1 


„ _3_ 1 „ _3_ S, 

iu m ^ { t*Sz 


^5= 3s 


3 S, 


+ n 


3 S^ __ 


« - 
3k 


S 3s 


^? 5 


_3 S,, 


>* 3k 
3 B, 


y* 


3y p. 


v J.s^ + V AM*. - v AE^_ 3_ B 
/®3k^ k1 >y 3y 'iy 3s 
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The proposed software suite involves the development of a systematic way of 1) computing 
derivatives, 2) carrying out TDMA, pentadiagonal solution (PDMA), or Jacobi/Gauss-Seidel 
iterations, and 3) updating the residuals. The schematic shown below represents the initial 
concept for the multidisciplinary software, while Appendix A shows a detailed design. The 
following remarks are in order. Note that within the “engine,” each application (CEM, MHD, 

FSI, CAA, etc.) could have, as an alternative to the “tt” formulations above, a separate routine 
that accesses the same derivative calculation function. The Phase I project focused on the “tt” 
formulation. An alternative, optimized approach will be explored in Phase II in order to obtain the 
most efficient code. This is intended for the applications that have less number of variables than 
the 8 in the base code. The graphics user interface (GUI) provides the appropriate software 
branches, as is the case in the current version of the commercial code INSTED. We will also have 
core and application-specific results and visualization (post-processing). 

Programming in Phase I was based on FORTRAN 90, C, and C++. The preliminary compute 
engine used FORTRAN, whereas C and C++ were used mostly for the interface codes. 

IV: Graphical User Interface (GUI) 

Although not promised for Phase I, we implemented a preliminary graphics user interface (GUI), 
to demonstrate the feasibility of a user-friendly interface for the proposed multidisciplinary 
software. This capability is an absolute necessity for commercialization of the final product. The 
GUI contains most aspects of the proposed software, except the engine (solver). That is, the left- 
hand side of the schematic shown below. Sample screens applicable to one of the applications are 
shown in Appendix B. 
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Figure 1 Preliminary Multidisciplinary Software Concept 
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V. Design of the multidisciplinary Software Suite 

A detailed design of the multidisciplinary software is shown in Appendix A. 

VI. Demonstration of the multidisciplinary Software Suite 


The feasibility of the proposed multidisciplinary software suite was established in Phase I, by 
ascertaining accuracy and computational efficiency when the proposed design in Appendix A is 
implemented. The base code for the development of the multidisciplinary code is the Air Force 
code FDL3D/MGD, which has a default number of variables of 8. We have modified this code to 
execute in a multidisciplinary environment and have used two problems for evaluation: 




MGD (Alfven wave with ohmic damping). Here, we compare the computational efficiency 
and accuracy of the original FDL3D/MGD code with those of the preliminary 
multidisciplinary code derived from it. The CPU times are shown in Table I and the accuracy 
below in figures 2 through 4. From the table we see that the efficiency of the 
multidisciplinary code is comparable to that of the original MGD code for both the Runge- 
Kutta and Beam-Warming time integration methods. Note that the CPU times shown are for 
100 time steps and the Beam-Warming calculations use 3 sub-iterations per time step. This 
exercise shows that the multidisciplinary software suite is quite feasible, in relation to the 
base code. 

Acoustic scattering by a concept X24C vehicle. This problem is interesting for a couple of 
reasons: it represents a complex geometry and involves only 5 variables (as opposed to the 8 
in the base code for the multidisciplinary software). Therefore, the associated overhead is 
expected to challenge the multidisciplinary code relative to a dedicated code for acoustic 
scattering. (See Ladeinde et al., 2001.) The initial results confirm our suspicion (not shown). 
However, a preliminary optimization technique was found to show great promise for the 
multidisciplinary software suite when the number of variables is much fewer than that in the 
base code. The latter will be fully developed in Phase II. 
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• Modification on Highlighted Blocks 



Modifications: 

1. In “Gcommon.f90” add code segment: 

i 

i 

i 

i 

i 

INTEGER:: MULTIDSCP 

INTEGER, DIMENSION^):: MULTI TYPE 

I 

I 

I 

I 

I 

2. In “I 


MULTI_TYPE=0 
READ *, MULTI DSCP 
MULTI_DSCP=MULTI_DSCP/10 
MULTI TYPE (MULTI_DSCP)= 1 
WRITE *, “WARNING SIGNS” 
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• Modifications on Transformation between Primitives and Conservatives 


PV to CV 



CALL CONVERT 
Or CALL DGCONVERT 

Or CALL UPDTQO 


CV to PV 



f 


L 


CALL UDGCONVERT 


COMMON CONTROL LOGICS: 


VARIABLEARRAY=0 ! INITIALIZATION 

| 

! CASE-SELECT STRUCTURE 
! 

MULITI ONE: SELECT CASE (MULTIDSCP) 

CASE (1:6) 

VARI ABLEARRA Y=BASIC TERMS in NS 
CASE (7) 

VARI ABLEARRAY=CEM TERMS 
END SELECT MULTI ONE 

MULITI TWO: SELECT CASE (MULTI DSCP) 

CASE (3) 

VARIABLE ARRAY=REACTING SPECIES 
CASE (4:5) 

VARI ABLEARRA Y=B FIELDS 
END SELECT MULTI TWO 


Notes: 

1. Replace the VARIABLE ARRAY by different WORK ARRAYS in different subroutines, 
namely, 

a) TEMP () in CONVERT.F90; 

b) RHS () in DGCONVERT.F90 

c) QOLD () in UPDTQO.F90 
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• Modifications on Beam-Warming Solver 
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CALL SWEEPJ1 

! (CALL RHSXZ) 

CALL SWEEPI 

! BY IVMON 

CALL SWEEPK 

! (CALL RHSY) 

CALL SWEEPJ2 
(CALL SDAMP/X/Y/Z/ 

in SWEEPK) 


CALL XICOMP 
COMPACT: CALL ETCOMP 

CALL ZTCOMP 


dF dG dH 


ROE: 


UPWIND: 


CALL XIROE 
CALL ETROE 

CALL ZTROE 

CALL XIUPWIND 
CALL ETUPWIND 

CALL ZTUPWIND 


dFy dG v dH v 


’ dr/ ’ 0 ^ 


COMPACT: 

CALL VCMPYRHS 
CALL VCMPXZRHS 

2 nd ORDER: 

CALL V2NDXZRHS 
CALL V2NDYRHS 
CALL VCROSX 


CALL VCROSY 

CALL VCROSZ 
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CALL CMATA 
CALL CMATB 
CALL CMATC 


5. 




CALL CMATV1 

a F v dG v dH v 


CALL CMATV2 

•) _ } 

<-> 

CALL CMATV3 

dQf dQ n dQ ( 




6. The same CASE-SELECT logics as used above are applied here. The work arrays specifically are 

a) WRK (*,*,*,*1), TEMP (*,*,*) and RHS (*,*,*) in SWEEP/I/J1/J2/K/.F90 

b) WRK (*,*,*,1) in /XI/ET/ZT/COMP.F90 

c) WRK (*,*,*,1) and RHS (*,*,*) in VCMP/XZ/Y/RHS.F90 and VCROS/X/Y/Z/.F90 

d) ALCL (*,*,*) in CMAT/A/B/C/.F90 and CMAT/V 1/V2/V3/.F90 

e) WRK (*,*,*,*) in SDAMP/X/Y/Z/.F90 

7. Note: Roe, Upwind, 2 nd -Order differencing schemes are currently not implemented for the multi¬ 
disciplinary computation. 

8. ADD SUBROUTINE MULTISOURCE.F90 for Beam-Warming procedure. The algorithm is: 


LHS of BW (0) 

At 

=RHS of BW (0) +- S 

1 + 0 

Where BW (0) represents the original Beam-Warming procedure and the 
expression of the source term S is described in my notebook. 
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• Modifications on Runge-Kutta Solver 



Detail Descriptions of the implementation: 


1 . 


Pre-Processing 


Post-Processing 


2 . 


CALL RUNGE 


25 





































DOD/SBIR Phase I Final Report 

Advanced Software Suite for Multidisciplinary Computation (Thaerocomp Technical Corp.) 


March 3, 2003 


r 



3. The same CASE-SELECT logics as used above are applied here. The work arrays, specifically, are 

a) TEMP (*,*,*) and QOLD (*,*,*,*) in RUNGE.F90 

b) WRK (*,*,*,1) in /XI/ET/ZT/COMP.F90 

c) WRK (*,*,*,1) and RHS (*,*,*) in VCMP/XZ/Y/RHS.F90 and VCROS/X/Y/Z/.F90 

4. ADD SUBROUTINE MULTISOURCE.F90 for Runge-Kutta procedure. The algorithm is: 


LHS ofRK(O) 

=RHS of RK (0) + A t -J S 

Where RK (0) represents the original Runge-Kutta procedure and the 
expression of the source term S is described in my notebook. 
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Appendix B: Sample GUI Screens 
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Table I: Comparison of CPU time (seconds) for the calculation of Alfven wave with ohmic damping. The 
CPU times are reported for 100 time steps. The Beam-Warming scheme uses 3 sub-iterations per time step. 


Code Type 

Grid Size 

CPU time (seconds) 



Runge-Kutta 

Beam-W arming 

Base code 

101x11x11 

3035.9 

4742.1 

Multidisciplinary code 

101x11x11 

2990.9 

4787.0 

Multidisciplinary code 

101x3x3 

3632.0 

95.09 
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Appendix C: List of FDL3D Files Modified 

The following is a list of FDL3d/MGD project files that were modified in the preliminary 
conversion of the code to a multidisciplinary suite. Note that some files contain multiple 
routines/procedures. 

bndry.f90, cmata.f90, cmatb.f90, cmatc.f90, cmatvl.f90, cmatv2.f90, cmatv3.f90, convrt.f90, 
datin.f90, deriv.f90, dgconvrt.f90, etcomp.f90, fdl3d_ice.f90, fdtry.f90, gcomonl.f90, 
met2nd.f90, multi_src.f90, nirhs.f90, output.f90, pirhs.f90, rhsrange.f90, range.f90, 
sdampx.f90, sdampy.f90, sdampz.f90, sdampzl.f90, sdiagx.f90, sdiagy.f90, sdiagz.f90, 
sweepi.f90, sweepjl.f90, sweepj2.f90, sweepk.f90, tmstep.f90, txiirhs.f90, tzetrhs.f90, 
udgconvrt.f90, updtqo.f90, vcmpxzrhs.f90, vcmpyrhs.f90, writetec.f90, wrtfield.f90, 
wrthist.f90, xicomp.f90, ztcomp.f90 
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APPENDIX D: Grid for Second Test Problem for the Multidisciplinary Software 




Fig. 5 Coarse grid model of X24C Reentry Vehicle showing the surface mesh. 
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